Self-similar correlation function in brain resting-state fMRI 
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Abstract 

Adaptive behavior, cognition and emotion are the resuh of a bewildering variety of brain spatiotem- 
poral activity patterns. An important problem in neuroscience is to understand the mechanism by 
which the human brain's 100 billion neurons and 100 trillion synapses manage to produce this large 
repertoire of cortical configurations in a flexible manner. In addition, it is recognized that temporal 
correlations across such configurations cannot be arbitrary, but they need to meet two conflicting 
demands: while diverse cortical areas should remain functionally segregated from each other, they 
must still perform as a collective, i.e., they are functionally integrated. Here, we investigate these 
large-scale dynamical properties by inspecting the character of the spatiotemporal correlations of 
brain resting-state activity. In physical systems, these correlations in space and time are captured 
by measuring the correlation coefficient between a signal recorded at two different points in space at 
two different times. We show that this two-point correlation function extracted from resting-state 
fMRI data exhibits self- similarity in space and time. In space, self- similarity is revealed by consid- 
ering three successive spatial coarse-graining steps while in time it is revealed by the 1// frequency 
behavior of the power spectrum. The uncovered dynamical self-similarity implies that the brain is 
spontaneously at a continuously changing (in space and time) intermediate state between two ex- 
tremes, one of excessive cortical integration and the other of complete segregation. This dynamical 
property may be seen as an important marker of brain well-being both in health and disease. 
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It is increasingly evident that brain regions are continu- 
ously interacting even when the brain is "at rest" and, 
more importantly, that the functional networks uncov- 
ered from resting data closely matches those derived from 
a wide variety of different activation conditions [H, 
Starting with the uncovering of coherent fluctuations of 
functional magnetic resonance imaging (fMRI) in time 
series of motor cortex Q, many other findings have val- 
idated the notion of correlated networks as a dynamical 
substrate of the resting brain. It has been established 
that these networks, which can be separated on the basis 
of their temporal features 0-0] , are located at consistent 
locations across subjects and are equally detectable even 
during sleep and anesthesia 



These exciting findings provide a novel window to ob- 
serve the brain at work and, at the same time, high- 
light our limited understanding of the functional orga- 
nization of the brain at large scales compared with 
the, often precise, knowledge we have of the small (neu- 
ral circuit level) scale. In particular, little is known on 
how the cortex is able to solve the conflicting dynam- 
ical demands imposed by the functional segregation of 
local areas differing in their anatomy and physiology on 
one side and on the other their global integration shown 
during perception and behavior. This riddle is clearly 
pointed out by Tononi fld\ : "traditionally, localizationist 
and holist views of brain function have exclusively em- 
phasized evidence for either functional segregation or for 
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functional integration among components of the nervous 
system. Neither of these views alone adequately accounts 
for the multiple levels at which interactions occur during 
brain activity" . In connection with this dilemma, it has 
been suggested that the brain's conflicting demands are 
a generic property of many collectives, regardless of be- 
ing composed by neurons, genes, individuals, etc [TTI.[l^. 
This implies that the search for the mechanism behind 
this dilemma must be guided by the general properties 
of the system, rather than by the details of the neurobi- 
ology. In that regard, the study of large scale collective 
properties have a long tradition in statistical physics, al- 
lowing the identification of different dynamical regimes 
through the study of ubiquitous correlation properties. 
In this work, we are guided along that direction, and 
the paper is dedicated to report the spatiotemporal cor- 
relation properties of fMRI resting-state data which are 
found to exhibit robust self-similarity signatures in space 
and time. This finding has profound significance, be- 
cause it demonstrates a continuous range of correlations 
between the local cortical circuits up to the entire brain. 
Thus, this is the first direct empirical demonstration of 
the brain resting-state activity exhibiting simultaneously 
functional segregation and integration. 

A variety of experiments have already provided indi- 
cations of self-similarity in spatial and temporal scales 
of the brain dynamics. For example, at small scale, 
avalanches of neuronal activity in rat [isl . IT^ and mon- 
key [IHl cortex are known to be scale-free, a finding that 
has been modeled by Arcangelis et al. Levina et al. 
[l7| and Buice and Cowan [3| among others. At larger 
scale, it was reported that brain fMRI networks are char- 
acterized by scale invariant degree distributions as 
described by Eguiluz and collaborators and later repli- 
cated and extended by van den Heuvel et al. [131 • These 
networks are indistinguishable from those extracted from 
model systems at criticality [H, [2lj . Other findings in- 
clude the observation of 1// power spectra from simul- 
taneously recorded magnetoencephalography and elec- 
troencephalography signals 22] , fMRI signals (23l - [25| as 
well as from cognitive responses [2(|] , all indicative of long 
range correlations in brain dynamics processes. It is these 
myriads of observations of power-law behavior and long 
range correlations that have led to the conjecture that 
the human brain as a whole behaves as a system at crit- 
icality f27|-p9l and that the framework of self-organized 
criticality HHllI may be re evant to understand large- 
scale brain dynamics. In that sense, critical dynamics 
endows the system with a high susceptibility and a broad 
repertoire of responses, natural requirements for healthy 
brain function. 

Here, we analyze directly the correlation function of 
the voxel signals in space and time. A voxel is assigned 
a set of three weights describing the density of grey mat- 
ter, white matter and cerebrospinal fluid (CSF) within 
the voxel. Up to now, the standard analysis of fMRI 
functional correlations only included voxels with a den- 
sity of grey matter above a certain threshold value. In the 



present analysis, we do not discard any voxels, nor do we 
apply a threshold on the correlation coefficients. Hence, 
the present analysis corresponds to the standard anal- 
ysis tool applied in physics and materials science when 
structural properties are investigated. 
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FIG. 1: The renormalized average correlation function 
2~""(C) vs real distance (inset: (C) vs voxel distance) for 
the four levels of description for one subject. Full line: 
128 X 128 X 31 (n = 0), dashed line: 64 x 64 x 16 (n = 1), 
dashed-dotted line: 32 x 32 x 8 (n = 2) and dotted line: 
16 x 16 x 4 (n = 3) . Top panel linear-linear and bottom panel 
log-log axis, respectively. Full straight line is a guide to the 
eye for a power law with exponent /3 = 0.47 ± 0.2. Parameter 
a = 0.45. 



To investigate whether the system is self-similar in 
space, we perform a coarse-graining analysis directly on 
the voxel data and then extract the associated correlation 
function. We will show that the coarse-grained correla- 
tion function exhibits a power-law decay of correlations 
at length scales between the microscopic scale of a voxel 
and the macroscopic scale of the size of the brain. This 
self-similarity is reminiscent of the multi-scale modular 
organization [H, observed in fMRI data [sl] , but also 
a robust finding in systems near the critical point [29l - [3l| . 



I. RESULTS 

Spatial Self- similarity: Starting from the computation 
of the correlation matrix, the average correlation func- 
tion for each separation r can be estimated. To study the 
degree of self-similarity present in the correlation func- 
tion, we perform a coarse-graining procedure and check 
for self-similarity via Eq. ([7]) as detailed in the section of 
Material and Methods. The result of this analysis for a 
single subject is shown in Fig. ([I}. The Supplementary 
Material (see below) contains the results for six addi- 
tional subjects. In all cases, the data were binned into 
integer distance values to decrease the statistical fluctu- 
ations (especially at long distances) without having to 
average across several subjects. Two different regimes 
are observed in the correlation function. At short dis- 
tances, we see a collapse of the different curves so Eq. ([7]) 
holds indicating that the system is self-similar. At longer 
distances, instead of the usual decay of the correlation 
function that is observed in physical models, there is a 
plateau, probably due to the brain bilateral symmetry, 
as seen previously across cortical bilaterally homologous 
regions by Salvador et al. [s^ [s^ , and then an increase 
in the correlations, probably due to surface effects. 

After one coarse-graining step, the length of a side of 
a block-voxel is twice the size of a voxel, but then we 
need to renormalize the Euclidean distance in the coarse- 
grained system by a factor 2 to recover distance in vox- 
els. If the system under consideration was infinite, then 
(C(r)) would be defined for all r irrespectively of the level 
of description and they would be identical, see Eq. ([T]). 
However, as the system we consider is finite, there exists 
a maximum distance in the original system. After 
one coarse-graining step, this maximal distance equal to 
r,„/2 and after n coarse-graining steps, this maximum 
distance equals rm/2", measured in units of voxels at 
the nth coarse-graining level. If Eq. ([7]) holds, then it 
means that the system behaves in the same manner at 
all scales, up to the maximum distance possible. Fur- 
thermore, it allows for a data-collapse of the correlation 
functions by first re the distance r in voxel units to Eu- 
clidean distances, r i— > r2"a where a is the voxel spacing 
in Euclidean distance (say mm) in the original data and 
then renormalizing the correlation functions by the fac- 
tor 2~"", where a is a subject-dependent parameter. It 
is easy to check that if the correlation function depends 
on distance through a power law with exponent /3, that 
is, (C(r)) oc r~^, then the exponent a in the renormal- 
izing factor of the correlation function must be equal to 
p. The relation a = /3 is fulfilled for this regime. 

Scaling is only expected to apply in an intermediate 
regime where the sizes of the coarse-grained voxels are 
big compared with the smallest length scale cut-off in the 
data, that is, the voxel scale set by the scanner, and the 
largest distance covered by the data, that is, the size of 
the brain. This is indeed observed in Fig. ([T]) where the 
data collapse is best for the three curves corresponding 
to the bigger block voxels. 
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FIG. 2: The average correlation function (C) vs the voxel 
distance after one step of coarse graining (64 x 64 x 16) for 
fMRI data (solid line) spatially randomized voxels (dashed- 
dotted line) and temporal randomized voxels (dotted line). 
Shown is also the average correlation function for phantom 
data (dashed line). Top panel in linear-linear, middle in 
linear-log and bottom in log-log axis, respectively. The cor- 
relation function for the fMRI data decays as a power law, 
while the correlation function for the phantom data decays 
exponentially fast. The correlation function for the spatially 
randomized data is constant while it is zero for the temporally 
randomized data. 



To verify that the shape of the average two-point cor- 
relation function is not due to artifacts in the data, we 
randomized the positions of the voxels or the time series, 
see Fig. ([2]). By randomizing the positions, we expect 
to find a constant profile for the correlations. In an infi- 
nite system or in a system where the correlation length is 
much smaller than the system's size, this constant would 
be equal to zero, however, in the system we study strong 
long-range correlations are present, hence the constant is 



non-zero. By randomizing the time series of each voxel 
independently, we expect to destroy all correlation among 
voxels and thus expect close to zero average correlation 
for all distances. Furthermore, we also scanned a phan- 
tom to test our finding against artifact coming from the 
scanner. The phantom data shows an exponential de- 
cay and no long range correlation, except a peak at the 
longest distance that is also observed in brain data, con- 
firming that it is a surface effect. 

Pondering the Correlation Function with the Voxels' 
Grey and White Matter Content: When analyzing fMRI 
data, the focus is, unlike our approach, generally on grey 
matter voxels. Due to the large millimeter size of the 
voxels, the discrimination between grey and white matter 
voxels is based on a high resolution anatomical scan that 
gives the content of grey matter, white matter and CSF 
for each voxels. A voxel is said to belong to grey/white 
matter if its content in the latter exceeds a threshold that 
has to be fixed, the others are dismissed. This procedure 
has the drawback that information is lost in the process, 
because the grey matter signal in dismissed voxels is lost. 

Pure grey and white matter voxels are rare, but from 
a high resolution scan, we know the content of grey and 
white matter of each voxel. Hence, we can weight the 
correlations between two voxel according to their rela- 
tive content of grey and white matter with the following 
normalization: 1 = + -\- w'~'^^ . We have per- 
formed the same analysis as we did above. Using all 
voxels weighted with their contents of grey and white 
matter, respectively, we calculated the average correla- 
tion functions 



(C-) = «(C,.,). (1) 



The correlation of grey matter is shown in Fig. |3] and 
the correlation of white matter is shown in Fig. 21 The 
bump at the end of the distribution survives only in the 
grey matter, which is to be expected, since it is proba- 
bly due to a surface effect. White matter correlations die 
off smoothly before the surface of the brain is reached, 
which was expected since white matter is supposed to 
be enclosed in grey matter. The fact that the exponents 
a 7^ /3 in both Figs. [3] and |3] is related to the deviation 
from power-law behavior of the correlation functions of 
the grey and white matter. The scaling collapse is better 
for the correlation function when all voxels are included 
(see Fig. [T]) than for the grey matter and white matter 
correlation functions displayed in Figs. [3] and [H It is 
remarkable that only when the entire fMRI BOLD signal 
(all voxels without thresholding) is included in the corre- 
lation matrix we do obtain a good quality data collapse 
consistent with self-similarity. In general, we associate 
the essential dynamics and processing of the brain with 
grey matter, so our findings indicate that the separation 
of the fMRI BOLD signal using the weights and 
is more subtle than expected. 

Temporal Self- similarity: To complement the spatial 
correlation analysis, we present the average power spec- 
tra of the fMRI time series. Temporal auto-correlations 
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FIG. 3: The renormalized average correlation function 
2~""(C®®) of signal weighted with the content of grey matter 
vs real distance (inset: (C"®) vs voxel distance) for the four 
levels of description for one subject. Full line: 128 x 128 x 31 
(n = 0), dashed line: 64 x 64 x 16 (n = 1), dashed-dotted 
line: 32 x 32 x 8 (n = 2) and dotted line: 16 x 16 x 4 (n = 3). 
Top panel linear-linear and bottom panel log-log axis, respec- 
tively. Full straight line is a guide to the eye for a power law 
of exponent /3 = 0.60 ± 0.2. Parameter a = 0.35. 



and power spectra have been studied by use of fMRI for 
a while. A comprehensive review is given in [37j and the 
use of these as a diagnostic tool was discussed in 2J] . To 



obtain the spectra, we computed the spectrum of each 
voxel individually and normalized its integral to 1. This 
corresponds to analyzing the Fourier transform of the 
auto- correlation function. We note that the term "power 
spectrum" is often used to denote the Fourier transform 
of the auto-covariance function. However, it is better 
to operate with the normalized correlation function to 
avoid signals with more power to totally dominate the 
average power spectra. Although the frequency range we 
can explore is rather small, the power spectrum has a 
reasonable 1// dependence, thus showing criticality in 
time as well as in space, see Fig. ([S]). To check against 
potential artifacts of the measurement, we calculated the 
power spectra of phantom data and surrogate data, con- 
structed by randomizing the brain fMRI time series. In 
both cases, we obtain a flat spectra as shown in Fig. (O 
rejecting the possibility that the spectral features arise 
from scanning artifacts. The peak at 0.03Hz superim- 
posed on the 1// slope corresponds to the low-frequency 
fluctuations dynamics described already for the so-called 
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FIG. 4: The renormalized average correlation function 
2~""(C"™) of signal weighted with the content of white mat- 
ter vs real distance (inset: (C™*') vs voxel distance) for 
the four levels of description for one subject. Full line: 
128 X 128 X 31 (n = 0), dashed line: 64 x 64 x 16 (n = 1), 
dashed-dotted line: 32 x 32 x 8 (n = 2) and dotted line: 
16 x 16 x 4 (n = 3). Top panel in linear- linear and bottom in 
log-log axis, respectively. Full straight line is a guide to the 
eye for a power law of exponent /3 — 0.61 ± 0.2. Parameter 
a = 0.20. 



brain resting-state networks Q. 



II. SUMMARY AND DISCUSSION 

These results are the first direct demonstration of spa- 
tial and temporal self-similarity of the brain resting state 
dynamics. Although several reports already hinted about 
the possibility that at large-scale the brain may be op- 
erating in such state, previous measures have been indi- 
rect and/or model dependent. In contrast, in the present 
work scale invariance in space and time is determined 
from the entire data set without discarding any voxels 
and using the same renormalization techniques champi- 
oned in the study of critical phenomena in physical sys- 
tems and models. 

It is interesting to relate our finding of scale invariance 
with the spatial and temporal localized functional activ- 
ity extracted from fMRI data by Independent Compo- 
nent Analysis (ICA) methods. By construction, ICA em- 
phasize the identification of strongly independent com- 
ponents and "filtering" out less important contributions 
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FIG. 5: Average power spectrum for the fMRI time series 
(solid line), the phantom time series (dashed line) and the 
randomized time series (dotted line). Full straight line is a 
guide to the eye for a 1// spectrum. 



to the total BOLD signal. Recent studies have investi- 
gated the growing complexity of the spatial structures 
extracted by ICA as the total numbers of components is 
allowed to increase 0, llHl • The additional structure ob- 
tained this way is indeed consistent with the self-similar 
character we have identified in the correlation function 
of the BOLD signal. 

It should be mentioned that the distribution of BOLD 
fMRI spatiotemporal correlations was found already to 
be altered in certain chronic diseases [39j which together 
with the present result seem to suggest that scale invari- 
ance could be an important objective marker of brain 
well-being both in health and disease. 

From a dynamical systems perspective, the uncovered 
self-similarity implies that the brain dynamics is perma- 
nently at an intermediate state between two extremes, 
one that is strongly correlated across large distances, pro- 
ducing transient highly integrated cortex states and the 
other in which only nearby clusters are acting in sync. 
This scenario, of long range correlations in space and 
time, is only conceivable in dynamical systems at criti- 
cality and could be the manner in which the cortex can 
manage to produce an arbitrarily large repertoire of inter- 
action patterns among arbitrarily distant cortical sites. 



III. MATERIALS AND METHODS 

Image Acquisition: A 3T Philips MRI scanner was 
used to acquire Tj-weighted echo-planer images (EPI) 
in 128 X 128 x 31 voxels of dimension 2.5 x 2.5 x 5.0 mm"^ 
using a repetition time of 2000 ms, echo time of 60 ms 
and a flip angle 90°. An 8-channel array coil and SENSE 
factor 2 were used as well as second-order shims. Finally, 
a 3-dimensional mask was used to identify the content of 
white matter, gray matter and CSF. 

Image Pre-processing: The brain data consist of fMRI 
data sets of seven young adults healthy subjects in the 
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resting state. The subjects were instructed to lie still in 
the scanner with their eyes closed avoiding falling asleep. 
A total of 305 functional volumes of each subject were ac- 
quired from each session, the first five of which were dis- 
carded in order to remove the effect of Ti equilibration. 
Image pre-processing was performed using the Univer- 
sity of Oxford's FMRIB software library (FSL) involving 
realigning to account for movement using FMRIB 's Lin- 
ear Image Registration Tool (MCFLIRT) [IQl, high-pass 
temporal filtering using FEAT to remove low frequency 
artifacts and a high pass filter cutoff preset to 100 s. 

Coarse Graining of the Signal: Coarse graining is a 
well established technique to describe system's changes at 
different levels of spatial or temporal observation. This 
allows one to investigate to what extent the spatial or 
temporal structures of a given phenomena exhibits scale 
invariance, see for example Ref. [ij]. A self-similar 
or scale invariance object, for instance a cauliflower, 
looks like itself at all scales. Each little bouquet of 
the cauliflower looks a miniature version of the entire 
cauliflower. Hence, a short length scale study of the 
cauliflower, in which one probes the substructure of an 
individual bouquet, will provide the same information as 
the study at a large length scale where one observe the 
cauliflower as a composition of bouquets. 

In general, the approach to probe the behavior at dif- 
ferent length scales consists in replacing the signal at a 
given point by a spatial average over a region of a given 
spatial extend about this point (see Fig. In the case 
of an fMRI signal, the procedure consists in aggregating 
2x2x2 = 8 adjacent voxels to create a block-voxel B' 
whose BOLD intensity at time t, VB'{t), is the average 
intensity of the BOLD intensities at time t, wb(0 in its 
constituting voxels B £ B' ^ that is. 



(2) 



see Fig. |6]for an illustration. 
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FIG. 6: (Left) Example of coarse graining in 2 dimensions 
where there are 4 boxes B within a block-box B' . (Right) The 
four dashed-colored signals from the four original boxes B are 
averaged to produce the solid-black coarse-grained signal of 
B'. 

Computing Eq. [5] for all times t = 1,2,..., 300 gives 
the time series for the BOLD signal in the block-voxels 



{visi{t)}. We now repeat the coarse-graining procedure 
on the time series for the BOLD signal in the block- voxels 
{vB'{t)} to produce a time series for the BOLD signal 
in the block- voxels {vt}"{t)} and so on. For the spatial 
resolution of our fMRI data, the coarse-graining proce- 
dure applied three times yield voxels ;B(1 28x1 28x31) H> 
e'(64x 64x16) ^ ^"(32x32x8) ^ B"'(16xl6x4) with 
associated signals VB{t) VBi{t) M- VB"{t) VB"'{t). In 
the first coarse-graining step, the ultimate plane of voxels 
at the 31*'' slice is coarse grained only along the x — y 
plane. Averages has been calculated over active voxels 
only to prevent the inclusion of non-active voxels at the 
boundaries. 

Next in our analysis, we consider the time series of each 
(block) voxel as a realization of a random variable and 
compute the time correlation matrix at the four levels of 
coarse graining. This matrix is obtained by computing 
the correlation between all pairs of (block) voxels Bi at 
position Ti and Bj at position rj : 



Cb {vb, , VBj ) = 



t=i 



(3) 



where /tg. — y X]t=i (0 i^ the temporal average 
value of (block) voxel's Bi BOLD signal, and — 

\J ^ Y^=i ^Bi (^) ~ /^Si standard deviation of the 

BOLD signal of (block) voxel Bi and T = 300. 

Scaling: From the correlation matrix, the correlation 
function at equal time (ue^, ug^. ) is calculated as a 
function of the position in space of the voxels. We then 
compute the correlation function (Ce(^)) averaged over 
all voxels with r = jr^ — rj| measured in units of (block) 
voxels. This is done for each level of coarse graining 
to obtain (Ce(r)), (Ce'(r)), {CB"{r)) and (Cf3"'(r)) and 
the indices B' , B" and B'" indicate the level of coarse 
graining 

Let us discuss what to expect concerning the relation of 
the correlation function calculated at one coarse-graining 
level, say (Ce(r)), compared to the one at the next level, 
{Cb' {r)). To clarify the efi^ect of coarse graining, we relate 
the correlation function at this level of coarse graining 
to the covariance function at the previous level. For a 
homogeneous system, the covariance function is given by 

Covs(r) = {vBiO)vBir)) - (t.s(0))(fs(0)), (4) 

where (•) denotes spatial and temporal average. We con- 
sider, like in Eq. ([2]), the coarse graining over a block- 
voxel B' centered at r: 



VB'ir) = -r^ VB{r + a), 



(5) 



r+aeB' 



where the block-voxel B' contains \B\ voxels. We then 
have 

o CovB(r -I- a2 — ai) 
(Cb,(0)= ^'"^ . . ■ (6) 
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For a scale invariant (self-similar) system of infinite size, 
{Cig>{r)) will remain invariant if distance r is measured 
in units of coarse-graining boxes, that is, 

{CB{r)) = (C^'ir)) = {Cs"{r)) = {CB'"{r)). (7) 

Note that r is the distance expressed in the voxel unit 
corresponding to the given level of coarse graining. 
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FIG. 7: Main: The average correlation function (C) vs the 
block distance for three coarse-graining levels of an hypothet- 
ical object with a power-law form of the correlation function. 
The graphs collapse in voxel distance. Inset: Results from 
an hypothetical object having an exponential form of the cor- 
relation function where the curves do not collapse in block 
distance. 



The simplest case of scale invariance corresponds to a 
power-law behavior of (C(r)). Figure [7] illustrates the 
different behavior under rescaling of two 3-dimensional 
artificial systems, where we imposed the correlation to 
have either a power-law or an exponential behavior: 

{Cir)) = (8a) 
I -I- r 

(C(r)) ==exp(-r/5). (8b) 

Note that only the correlations having a power-law form 
do collapse in block unit. 
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Supplementary information 

The three average correlations (C), (C^^) and (C""") at four levels of description for six additional subjects are 
consistent with the those presented in the main text. In general we find that the equality q = /3 is more accurately 
satisfied when only the three most coarse-grained data sets are included in the data collapse. Full line: 128 x 128 x 31 
(n = 0), dashed hnc: 64 x 64 x 16 (n 1), dashcd-dottcd hnc: 32 x 32 x 8 (n = 2) and dotted line: 16 x 16 x 4 
(n = 3). Full straight line are guides to the eye for a power law with exponent /3. 




Distance in mm Distance in mm 





FIG. SI: The renormalized average correlation functions 2 ™"(C) vs real distance (inset: (C) vs voxel distance). Left panel 

linear-linear and right panel log-log axis, respectively. Top panel 2~""{C) with /3 = 0.64 ± 0.2 and a — 0.51, middle panel 
2-'""(C""') with 13 = 0.74 ± 0.2 and a = 0.50 and bottom panel 2-°'"{C""") with /3 = 0.65 ± 0.2 and a = 0.40, respectively. 
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FIG. S2: The renormalized average correlation functions 2 ""(C) vs real distance (inset: (C) vs voxel distance). Left panel 
linear- linear and right panel log-log axis, respectively. Top panel 2~""(C) with /3 = 0.51 ± 0.2 and a = 0.42, middle panel 
2^""(C""') with /3 = 0.62 ± 0.2 and a = 0.37 and bottom panel 2-°"(C"'"') with ^ = 0.65 ± 0.2 and a = 0.31, respectively. 
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FIG. S3: The renormalized average correlation functions 2 ™"(C) vs real distance (inset: (C) vs voxel distance). Left panel 
linear- linear and right panel log-log axis, respectively. Top panel 2~°'"(C) with /3 = 0.46 ± 0.2 and a = 0.37, middle panel 
2-""(C'"') with /3 = 0.55 ± 0.2 and a = 0.40 and bottom panel 2-«"(C"'"') with = 0.52 ± 0.2 and a = 0.25, respectively. 
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FIG. S4: The ronormalizcd average correlation functions 2 ""{C) vs real distance (inset: (C) vs voxel distance). Left panel 
linear-linear and right panel log-log ikxis, respectively. Top panel 2~""{C) with /3 = 0.50 ± 0.2 and a = 0.41, middle panel 
2-""(C»») with /3 = 0.60 ± 0.2 and a = 0.37 and bottom panel 2-«"(C"'"') with /3 = 0.57 ± 0.2 and a = 0.30, respectively. 
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FIG. S5: The renormalized average correlation functions 2 ""{C) vs real distance (inset: (C) vs voxel distance). Left panel 
linear-linear and right panel log-log axis, respectively. Top panel 2~""(C) with /3 = 0.47 ± 0.2 and a — 0.34, middle panel 
2-""{C^^) with P = 0.72 ± 0.2 and a = 0.58 and bottom panel 2-""{C""'") with /3 = 0.55 ± 0.2 and a = 0.17, respectively. 
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FIG. S6: The renormalized average correlation functions 2 ™"(C) vs real distance (inset: (C) vs voxel distance). Left panel 
linear- linear and right panel log-log axis, respectively. Top panel 2~°'"(C) with /3 = 0.49 ± 0.2 and a = 0.39, middle panel 
2-""(C'"') with /3 = 0.60 ± 0.2 and a = 0.36 and bottom panel 2-«"(C"'"') with = 0.62 ± 0.2 and a = 0.29, respectively. 



